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Abstract 

We study the limited data problem of the spherical Radon transform in two and three 
dimensional spaces with general acquisition surfaces. In such situations, it is known that the 
application of filtered-backprojection reconstruction formulas might generate added artifacts 
and degrade the quality of reconstructions. In this article, we explicitly analyze a family 
of such inversion formulas, depending on a smoothing function that vanishes to order k on 
the boundary of the acquisition surfaces. We show that the artifacts are k orders smoother 
than their generating singularity. Moreover, in two dimensional space, if the generating 
singularity is conormal satisfying a generic condition then the artifacts are even k + I orders 
smoother than the generating singularity. Our analysis for three dimensional space contains 
an important idea of lifting up a space. We also explore the theoretical findings in a series of 
numerical experiments. Our experiments show that a good choice of the smoothing function 
might lead to a significant improvement of reconstruction quality. 


1 Introduction 

Let S C M n be a convex closed smooth hyper-surface. We consider the following spherical Radon 
transform 7Zf of a function / defined in M n 

TZf(z,r)= j f(y)da(y), (z, r) € S x (0, oo). 

S(jt,r) 


Here, §(z, r) is the sphere centered at z of radius r, and da is its surface measure. This transform 
appears in several imaging modalities, such as thermo/photoacoustic tomography (e.g., FPR04 
IFHR071lKK08| h ultrasound imaging (e.g., [NL791 lNL81| h SONAR (e.g., |QRS11| ), and inverse 
elasticity (e.g., fBK78j l. For example, in thermo/photoacoustic tomography (TAT/PAT), / is 
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the initial ultrasound pressure generated by the thermo/photo-elastic effect. It contains useful 
information about the inner structure of the tissue, which can be used, e.g., for cancer detection. 
On the other hand, the knowledge of 1Z(f)(z,.) can be extracted from the ultrasound signals 
collected by a transducer located at z € <S, which is called the observation surface. One, 
therefore, can concentrate on finding / given IZ(f). The same problem also arises in other 
aforementioned image modalities. 

It is commonly assumed that / is supported inside the bounded domain 0 whose boundary 
is S. Let us discuss an inversion formula under this assumption. Let V : C'q°(M+) —>■ C' 00 (M_|_) 
be the pseudo-differential operator defined by 

V(h)(r) = J fe^\Xr'h(s)dsdX, (1) 

M R+ 


and B : C°°(S x M + ) —>• C°°(fl) be the back-projection type operator 
B(ff)(x) = 2 ~^ J (z-x,u z ) g(z,\x - z\)da(z). 

S 

When S is an (n — l)-dimensional ellipsoid, one has the following inversion formula |Kun071 
INat.121 lHalli] H 

f(x) = BVlZf(x ), for all x G 17. (2) 

We note here that formula m was written in other forms in the above references. The above 
form, presented in |Ngul3| , is convenient to analyze from the microlocal point of view. Another 
advantage of the above form is that it can be implemented straight forwardly: 1) V can be 
computed fast by using Fast Fourier Transform (FFT) and 2) B only involves a simple integration 
on S. Fig. |T] is the result of our implementation when n = 2 and S is the unit circle. The image 
size is N = 2048 pixels. The sampling data has the resolution of n a = n r = 2048 for the spatial 
(angular) variable 2 = (cos 9, sin 0) and radial variable r € [0,2], The reconstruction is almost 
perfect. 

When <S is a general convex surface, the operator on the right hand side of (J2J) might not be the 
identity. However, it onlv differs from the identity bv a conmact ODerator (see |Natl2, Ha ED, 
which is a pseudo-differential operator of order —1 (as shown in |Ngul3| ). Moreover, from the 
numerical experiments, we observe that BVlZf and / are very close, even when S is not any 
ellipse/ellipsoid. For example, in Fig. [2] we present the reconstruction using BV1Z for S being 
the polar curve defined by 

f i 1 

<(x,y):x = -((2 + cos 0) cos 9 — l), y = -(2 + cos 9) sin#), 0 < 9 < 27r 

I Z Zi 

The reconstruction is almost perfect. 

lr The reader is referred to, e.g., FPR041 |FHR07 ( i PaI12[fSa!12l IPall4l for other inversion formulas. 
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(b) Reconstruction 




(c) Line profile along the central 
horizontal line: the phantom (left) 
and reconstruction (right) 


Figure 1: Reconstruction using the operator BV1Z (for n = 2) when S is the unit circle. Here, 
the number of angular samples as well as the number of radial samples of the data 7 Zf was 
chosen as n a = n r = 2048. 





(c) Line profile along the central 
horizontal line: the phantom (left) 
and reconstruction (right) 


Figure 2: Reconstruction using the operator BTTZ (for n = 2) when S is a polar curve. Here, the 
number of angular samples as well as the number of radial samples of the data IZf was chosen 
as n a = n r = 2048. 


In this article, we are interested in the limited data problem. That is, the knowledge of 
IZf is only available on a closed proper subset T C S with smooth boundary 3T and nontrivial 
interior Int(F) (see, e.g., [XWAK091 |QRSlTj ). It is natural to consider the following formula 

T 0 f(x) := B X rVKf(x ), (3) 

where xr is the characteristic function of T. As we observe from Fig. [3l formula (|3|) reconstructs 
some singularities and also smoothens out some singularities of the original image. Moreover, 
it also creates some artifacts (added singularities) into the picture. More detailed discussion 
will be presented in Sections [2] and [3] (see also Section [4] for numerical demonstrations). The 
operator To reconstructs all the visible singularities. However, the artifacts it generates are quite 
strong. We now introduce a generalization of To in order to reduce the artifacts. Let us consider 

Tf{x) := B X VUf(x). (4) 
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(b) Reconstruction 



Figure 3: Reconstruction from limited view data, collected on a quarter of a unit circle, using 
the standard reconstruction operator 7o with no artifact reduction. The acquisition surface is 
illustrated by the green line in the phantom image. 


Here, x £ C 00 ^) and x = 0 on S \ T. Moreover, we assume further that x > 0 °n the interior 
Int(r) of r and x vanishes to order k on the boundary <9T. Of course, T = 7o if x = XT- 

In this article, we will study formula dU). Namely, we will analyze which singularities are 
reconstructed and how strong the reconstructed singularities are, compared to the original ones. 
Moreover, we will describe how the artifacts are generated by (HJ) and how strong they are. 

Let us discuss here the main idea of our analysis. We first follow the approach in | |Ngul3| to 
represent T as an oscillatory integral. To that end, let y, be the Schwartz kernel of T. Then, it 
can be written as (see |Ngul3| ) 

M®, y) = ^2 / / / l 2 ) A I Ar 1 (z - x, v z ) x{z) d\ da(z). (5) 

r r R n 


By the simple change of variables (z, A) —>• £ = 2(z — x) A, 


K x ,y) = 


2(2 tt) 


1 _ [ A< x -V'&+\ x -v\ 2 2\xJ%,t)\ ) 


+ 


/ 

R n 

wl 


J.{[x-y,£)- \x—y\ 


2(2v 

R n 

= li+(x,y) + n-(x,y). 


2 HI 

2|:r-z_(x,£)| 


x(z+)df 


x(z+) df 


( 6 ) 


Here, z± = z±(x, f) are the intersections of S with the positive and negative rays 


-R+(s,£) = {(x + tO ■ t > o}, 

R-(x,£) = {(x + tf) : t < 0}. 


To illuminate the idea, let us at the moment assume that x vanishes to infinite order at the 
boundary of T C S, then the functions a±(x,£) = x( z ±) are (smooth) symbols of order zero. 
Then, due to |Sog93 Theorem 3.2.1], T is a pseudo-differential operator (TDO) of order zero 
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whose principal symbol is (see |SU09] . and also |SU13] . for the same result in a more general 
framework): 

<ro(z,£) = ~[x(z + )+ x(z-)\. 

One important consequence of T being a TDO is that it does not generate the artifacts into the 
picture. Moreover, let us denote £(x,£) = R+(x,£) U i?_(x,£), the line passing through x along 
direction £. Then, the following discussions hold. 

a) Let us denote 


V = {(x,£) G T*Q \ 0 : £(x,£) intersects Int(T)}. 

Then, any singularities (x, £) of / in this zone generates a corresponding singularity of 
on the observed data g := TZ(f)\r- It is, therefore, called the visible zone (see, e.g., 
|LQ00| IPalOOl IXWAK09 ] ). Let (x,£) G V, then either x( z +) > 0 or x( z ~) > 0- Therefore, 
(j 0 (x, £) > 0. That is, T is an elliptic pseudo-differential operator of order 0 near (x, £). Due 
to the standard theory of pseudo-differential operators (see Lemma lA.l0D . (x,£) G WF S (/) 
if and only if (x,£) G WF S (T/). That is, all the visible singularities of / are reconstructed 
by T with the correct order. We notice that one visible singularity may be either visible 
in two directions, when both z+(x,£) and z-(x,£) belong to Int(T), or in one direction, 
when only one of z+(x,£) or z_(x,£) belongs to T. We will say that they are doubly 
visible and singly visible, respectively. 

b) On the other hand, let us denote 

X = {(x,£) G T*fl \ 0 : £(x,£) does not intersect T}. 

Then, X is called invisible zone, since any singularity of / in this zone does not generate 
any singularity of g. A singularity of / in this zone is called invisible. We notice that 
for each (x,£) G X, then x( z ±) = 0- That is, the full symbol of g is zero near (x,£). 
Therefore, due to the standard theory of pseudo-differential operators, (x,£) WF (T/), 
even if (x,£) G WF(/). That is, all the invisible singularities are completely smoothen out 
by T. 

Let us now consider the case of our interest: x on ly vanishes to a finite order k on the 
boundary of T (k can be zero as in case of 7o). Then, T is not a TDO anymore, since a± are not 
smooth. Therefore, T may generate artifacts into the picture, as shown Fig. [3] where T = To- 

Our results (see Sections [2] and [3]), show that the effect of T on the zones V and X is exactly 
what we discuss above for the case of infinitely smooth y. This can be seen from the facts 
that a±(x,£) are smooth on these two zones. The artifacts, on the other hand, come from the 
boundary zone where at least one of a+ or a_ is non-smooth: 

dV = {(x,£) G T*n \ 0 : £(x,£) intersects <9(r)}. 

We will characterize how these artifacts are generated and how strong they are. To that end, 
we will make use of the formula ([5]). Compared to (|6|), the formula flSJ) still keeps track of the 
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geometric information of 1Z, which is helpful to understand the generation artifacts. Moreover, 
we will study ([5]) as a Fourier integral operator, whose order determines the strength of artifacts 
(compared to the original singularity generating them). 


In term of geometric description of the artifacts, we will take advantage of the technique 
developed in |FQ15| (see also |FQ14| ). In |Ngul5b| , the third author studied the strength of the 
artifacts in S is a hyperplane (that is, T is flat). In this article, we study the problem for any 
convex smooth surface S. We will restrict ourselves to the two and three dimensional problems 
(i.e., n = 2,3), since they are the most practical ones. We will follow the microlocal analysis 
technique in |Ngul5b . However, due to the generality of the geometry involved, the arguments 
are more sophisticated. Moreover, for the three dimensional problem, we introduce a new idea 
of lifting up the space, which is simple but interesting. It is worth mentioning that similar 
problem has been studied for the X-ray (or classical Radon) transform [ Kat97l 1KR921IRK921 
FQ13[ [Ngul5a| . 


The article is organized as follows. In Section [2l we state and prove the main results for the 
two dimensional problem. In Section [3l we state and prove the results for three dimensional 
problem. We then present some numerical experiments for the two dimensional problem in 
Section [dj Finally, we recall some essential background in microlocally analysis in the Appendix. 


2 Two dimensional problems and some numerical demonstra¬ 
tions 

Let us consider n = 2. We assume that S is a smooth convex curve, and T is a connected piece 
of S. For our convenience, we assume that T is arc-length parametrized by the smooth function 

z = z(s ) : [a, b] —>• M 2 , 

and its end points are a = z(a ) and b = z(b). Let us define the function h = h{s) : [a, b\ —>• M by 

K*) = X(z{s)). 

Then, h(s) vanishes to order k at s = a and b, and h(s) > 0 for a < s < b. 

We define the following canonical relations in (T*0 \ 0) x (T*0 \ 0) 

Ay = {(x,£;x,£) : (x,£) € V}, 


and 

A a = {(x, t(x — a); y, r(y — a)) € T*0 x T*0 : \x — a| = \y — a|, r ^ 0}, 

Ab = {(x, t(x — b); y, r{y — b)) € T*0 x T*0 : \x — b| = \y — b|, r ^ 0}. 

We notice that (x, £; y, y) £ A a iff (y, rj) is in the boundary zone, corresponding to a, and (x, £) is 
obtained from (y, rj) by rotating around the corresponding boundary point a. Similar description 
holds for Ab- 

Proposition 2.1. We have 

WF{y)' C A v U A a U A b . 
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Here, WF(/j)' is the twisted wave front set of p: 


WF (p)' = {(x, £ y, —rj) : (*, y, v) € WF( M )}. 

Proposition 12.II was proved in [Ngul5b| when T is a line segment. The proof carries naturally to 
the general curve T without any major changes. We present it here for the sake of completeness 
and convenience in later discussion. 


Proof of Proposition 12. 1L Due to the composition rule for wave front sets (see Theorem IA.4I 
in Appendix), we obtain □ 

WF(//)' C WF (n B y o WFixvmY, (7) 

where pg and p-p are the Schwartz kernel of B and 1Z, respectively. Let us now proceed to 
analyze the right hand side of the above inclusion. 

We note that 1Z is an FIO with the phase function (see, e.g., |LQ00[ PallOl, [Ngul3| ) 

(j)(z, r, x , A) = (|x — z | 2 — r 2 ) A. 


For the sake of simplicity, for an each (z,r) = ( z(s),r ) E S x M + , we use the notation 
(z,r,p,q) for ( z,r,pds,q dr ) E T*(«S x M + ). Due to Theorem IA.11I (see Appendix), we obtain, 
by letting t = 2A 


WF(^) C C n := {{z = z(s),r, t (x - z(s), z'(s )) , -rr; x,t(x- z)) : 

s € M, r G M + , x E D, \x — z\ = r, r ^ 0}. 


Also considering x( z ) as a function of ( z,r,x ), we have 

WF(x) C {(a, r, r',0; x,0) : r / 0} U {(b,r, r',0; x, 0) : r / 0}. 

Applying the product rule for wave front sets (see Theorem IA. 31 in Appendix), we obtain 

WF (xBn) C Ctz, x u C a U C 6 , (8) 


where 

Cn,x = C n n {z E supp(x)}, 

C a = {(a, r, t — a, z'(a)} T t ', — r r; x, r (x — a)) : |x — a| = r, t 7 / 0}, 

Cb = { (b, r, r (x — b, z'(b)) + r', — r r; x, t {x — b)) : \x — b| = r, t' / 0}. 

On the other hand, we notice that H is a FIO with the same phase function cj) as 1Z (but the 
order of variables is switched), see e.g., |Ngul3| . Therefore, 

wf( mb ) c dk, 

2 Since V is a pseudo-differential operator, it does not generate new wave front set elements. This well-known 

pseudo-locality property of a pseudo-differential operator is recalled in Appendix lA.21 see (1261) . 
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where C l is the transpose relation of C 


C* = {{x,^,s,r,p,q) : (s,r,p,q-,x,£) G C}. 


From (ED , we arrive to 

WF(/i) C (Cft)' o (Ctz jX U C a U Cb) 7 = (C 7 ?.) ° (Ciz,x u C a U C b ) = (C\i ° Ctz jX ) U (C^> o C a ) U (C-^ o Cb). 
We notice that 


Therefore, 


Cyj o — Ay, Cji o C a — A a , Cji o Cb — Ab- 
WF(/i)' C Ay U A a UA b . 


□ 


Remark 2.2. Assume that x( z ) = 0 in a neighborhood of zq G Int(T) and (x, £) G T*9 such 
that £ = r(x — zo). Then, due to 

(®,0 £ tt r ( WF ( x ^)), 


where i tr is the right projection operator. From we obtain: 


(x,bx,£) 0 WF(^). 


This observation will be used later in the proof of Theorem \2.3\ a). 


Let us now employ Proposition 12.11 to describe the geometric effects of T on the wave front 
set of / (see also the discussion in |FQ15|). We first keep in mind the following inclusion, coming 
from Theorem IA.21 

WF(77) C WF(/x)' o WF(/). 


Therefore, due to Proposition 12.11 


WF(T/) C [Ay o WF(/)] U [A a o WF(/)] U [A b o WF(/)]. (9) 

The first part on the right hand side contains all the singularities that may be possibly recon¬ 
structed by T. The other two contain all the possible artifacts generated by T. We now discuss 
the implications of © in more details. 

a) Let (x,£) G WF(/) be an invisible singularity. We observe that 


Ay o (x, £) = 0, A a o (x,£) = 0, A b o(x,£) = 0. 


Therefore, from inclusion ([9]) , (x,£) is not reconstructed and does not generate any arti¬ 
facts. 
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b) Let (x,£) £ WF(/) be a visible singularity. Then, 

Ay 0 (x,0 = (x,0, Aa°(z,£) = 0, A b o(x,£) = 0. 

From inclusion ([9]), (x,£) may be reconstructed and does not generate any artifacts. 

c) Let (x, £) £ WF(/) be a boundary singularity pointing through a, that is £ = r{x — a) for 
some r ^ 0. Then, 

A a ° {x, 0 = {(V, V) ■ \x - a| = \y - a|, rj = r(y - a)}. 

From the inclusion (j9}|, we observe that (x,£) may generate artifacts (y,rj) by rotating 
around a. 

Conversely, assume that (y,rj = r(y — a)) be an artifact. Then, there is (x,£ = t(x — a)) £ 
WF(/) that generates (y,y). 

Similar description holds for a boundary singularity pointing through b. 

The strength of the reconstructed singularities, described in b), will be obtained by analyzing 
the singularities of y near A \ (A a U Ab). This will be done by the standard theory of pseudo¬ 
differential operators. In order to analyze the strength of the artifacts, described in c), we will 
make use of a class of FIOs associated to a point, introduced in Section lA.3.11 

Here is our main result of this section: 

Theorem 2.3. We have 

a) Microlocally on A \ (A a U Ab), we have y £ I°( A). Moreover, its principal symbol is 

Mx, Q = \ bcM + xO-)], 

where z± is the intersection of the ray R±(x, f ) with S. 

b) We can write y = y a + yb such that: 

i) WF(y a ) C A U A a . Moreover, microlocally near A a \ A, y a £ I~ k ~ 2 (A a ). 

ii) WF(yb) C A U Ab- Moreover, microlocally near Ab \ A, y^ £ I~ k ~ a(Ab). 

Will present the proof of Theorem 12.31 in Section T2.11 We now describe some of its consequences. 

1) Let (x,£) £ WF(/) be a visible singularity. Then, due to Theorem 12.31 bh microlocally 
near (x,£;x,£), y is a Fourier distribution of order zero with positive principal symbol. 
Applying Lemma lA. 101 we obtain (x,£) £ WF S (/) if and only if (x,£) £ WF s (Tf). That 
is, all the visible singularities are reconstructed with the correct order. 

Moreover, the formula oo(x,£) provides the magnitude of the main part of reconstructed 
singularities. For example, if (x,£) is a jump singularity across a curve S with the jump 
equal to 1, then (x,£) is also a jump singularity across S with the jump equal to <to(x,£). 
This explains the difference in the magnitude of the reconstructed singularities, that is 
demonstrated in Section |4j 
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2) Now, assume that (x,£) £ WF S (T f) is an artifact pointing through a. Then, each of its 
generating singularities (y,rj) £ WF(/) satisfies 

(x,&y,v) € A a- 

Then, due to Theorem 12.31 b) and Lemma lA. 151 at least one generating singularity (y,r/) 
satisfies ( y,r]) £ WF s _fc(/). That is, all the artifacts are at least k order(s) smoother then 
their strongest generating singularities. 

If we assume further that (x. £) has finitely many generating singularities ( y,r}), each of 
them is a conormal singularity of order r along a curve whose order of contact with the 
circle S(a, \y — a| = \x — a|) is exactly lEI- Then, due to Theorem 12.31 bl and Lemma lA.161 
(x, £) £ WF(T/) is a conormal singularity of order r + k + \ along the circle §(a, \x — a|). 
That is, the artifacts are at least (k + order(s) smoother than the strongest generating 
singularity. In our numerical experiments in Section 01 we will demonstrate this fact. 


2.1 Proof of Theorem 12.31 

Let us now discuss the proof of Theorem 12.31 It is similar to that of |Ngul5b , Theorem 2.2]. 
However, we need to employ more sophisticated microlocal arguments due to the generality 
of the geometry involved. As in |Ngul5b Theorem 2.2], the proofs for a) and b) require two 
different oscillatory integral representations for y. 


Proof of a). The idea is similar to the case of infinitely smooth x presented in |Ngul3| (see 
also jSU09] for more general framework). The main point here is to microlocalize the argument 
to stay away from An(A a nAb). Let (x*,£*) £ T*fl be such that (x*, £*; x*, £*) £ A\(A a UAb). 
We need to prove that there is fx o £ I~ k ~ 2 (A) such that 


(x*,f>*,r)0WFGu- M ) (10) 

and the principall symbol of /io at (a;*,£*) equals |[x(Af) + x(zL)], where z± = z± (x*, £,*). 

Since (x*, x*, £*) € A \ (A a U Ab), % is smooth near both z = z\ and 2 r* . Let us define 

a cut-off function c such that c(z) = 1 near z = z± and zero elsewhere such that c \' £ C°°(S). 
Let us write 

T:=Tc + T', (11) 

where 

T c := £> r c x VTZ, r ■■= B T (1 - c) X VU. 

Let y! be the Schwartz kernel of T'. Since \(1 — c) = 0 at z = z±, using Remark 12.21 we 
obtain 

(x*,e-,x*,n? wf(/). 

3 We note here that the contact order is always at least 1, since both curves are perpendicular to rj at y. 
Therefore, the condition on the contact order is quite generic. 
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On the other hand, similarly to ([ 6 |) (see also |Ngul3| for the derivation), we obtain the 
formula for the Schwartz kernel of T c 

Hc( x ,y) = 2 J e^ {x v ’® + ^ x yl ai«,-, +( x,£)i ) c(z+)x(z+) d£ 

R 2 


We notice that a±(x,£) := c(z±) x(z±) is smooth and homogenous of degree 0 in £. Therefore, 

y c g /°(A). 

Moreover, the principal symbol of y is (see, e.g., |Sog93 Theorem 3.2.1]) 

co,c(ab£) = c(z + )x(z + ) + c(z-)x(z-). 


Since c(z) = 1 at z = z±, we obtain 

^o,c(a:*,D = x{z*+) + x{z*-)- 
This finishes the proof (HOl) where /iq = y c . 

Proof of b). Let us decompose x into the form 


X = Xa + Xb, 

where % a , Xb G C°°(T) such that Xa vanishes near b and Xh vanishes near a. We then can write 
T = Ta + % where 

Ta := By XaVlZ, Th ■= Br XbVIZ. 

Repeating the argument in Proposition 12.11 (see also Remark 12. 2p . we obtain 

WF(/i a ) C A U A a , WFWCAU4. 

We now prove that microlocally near A \ A a , y a G I~ k ~ 2 (A a ). Let (x*, £*; y*, rf) G A a \ A. 
We need to prove that there is yo G I~ k ~ 2 (A a ) such that 

WF(y a -y 0 ). (12) 

Indeed, we notice that x* and y* belong to the same circle centered at a. Therefore, 

(x* -y*,z'(a)) ± 0. 

Otherwise, the tangent line of T at a passes through the midpoint of the line segment connecting 
x* and y*. This would be a contradiction to the assumption that hi is convex and x*,y* G fb 

Thus, there is a neighborhood U of (x*,y*) and e > 0 such that 

(x — y, z'(s)') 7 ^ 0, for all (x, y) G U, sG [a, a + e\. 
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Let c € C°°(S) such that c = 1 near z = a and c = 0 for z = z(s) such that |s — a| > e. Let 
us write 

Ta := Tc + T', (13) 

where 

Tc : = S r c X a VIZ, T' := By (1 - c) Xa VV. 

We denote by fi c and \j! the Schwartz kernels of T c and T', respectively. An argument similar 
to the proof of Proposition 12.11 shows that 

WF(/i c ) cAuA a , WF(ju') C A. 


Therefore, 

Similarly to ([5]) . we can write: 

1 


Me (x,y) = 


2 7T 2 




; i(|,-d 2 -P-d 2 )A | A | {z _ ^ Vz) c{z) Xa(z) dX d(j{z y 


r r i" 

The phase function of /r c can be written as 

(|x - z\ 2 - |y - z\ 2 )\ = (x - y, x + y - 2z) A. 

Therefore, 


Me (x,y) = 


1 


2 7T 2 


s i( X -y, X+ y-2z)X\ x \ ( z _ Xj]/g) c (z) Xa (z) d Z d\ 


J_ J e i{x~y,x+y)X | A | J e —i(x-y,2z)\ 


(z - x, v z ) c(z) Xa (z) dz d A. 


Let us denote 


A{x,y,\) = J 


= / e -^ x -y^ z ) x i? - 


(z-x,v z ) c(z) Xa (z)dz. 


Applying Lemma 12.41 we obtain 

Me (X,y) = _i_| e *[(^^)-<-^a)]A |A| 


e -i(x-y, 2 a)A A (x, y, A) 


d A 


_L J e i(\x-*\ 2 -\y-*\ 2 )x B {x,y,X)d\, 


where B{x, y , A) is a symbol of order —(k + 1) on U. We notice that A a is the canonical relation 
associated to the phase function of y c (see Appendix IA.31) . Therefore, 


Mel u e I k 2 (A a ). 


12 











Let ho be obtained from y c by multiplying with a cut-off function near (x*, £*; y*, rj*). We arrive 
to no £ I~ k ~ 5 (A a ) and 


£ WF(/x a -/x 0 ). 


We conclude that /x a E / 2 (A a ) microlocally near (x* ,f*;y* ,ri*). 

The proof for /Xb is similar. 

In the above proof, we have used the following result: 

Lemma 2.4. Let h £ C' 00 [a,&] and define 


□ 


A(x,y,X) = J e- i < x ~ y ’ 2z ^ x h(s)ds. 

a 

We have A £ C 00 ^ x fl x R). Moreover, 

a) Assume that h(s) = 0 near s = b and h vanishes to order k at t = a . Then, on the set 
{(re, y) £ fl x Ll : (x — y, z'(s )) 7^ 0, for all s £ [a, b] such that z(s ) € supp(/i)}, 
we have 

—i {x—y, 2 a)A .*/„ .. w h^ \a) 


A(x,y, A) 


[l + r(x,y,A)]. 


[2i (x — y,z'(a )) A]( fc+1 ) 
ifere, z/ a = v{a) is the normal vector of S at a = z(a). 
b) Assume that h(s) = 0 near s = a, and h vanishes at t = b to order k. On the set 

{(x, y, A) € Ll x 17 x R : (x — y, z'(r)} 0, for all t £ [a,b\, z(s) £ supp(/t)}, 

we have 

—h^ k \b) 


-i(x-y,2b) X A ( A1 ~ - 

1 ,y ’ j [2* (x-y,z'(b)) A]( fc +D 
Here, = 1 fib) is the normal vector of S at b = z(b ). 

In both a) and b), r(x,y,X) is a symbol of order at most —1. 


[l + r(x,y, A)]. 


Proof. The Lemma is proved by successive integration by parts. It is very similar to |Ngul5 b 
Lemma 2.3]. We skip it here for the sake of brevity. □ 
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3 Three dimensional problem 

Let us now consider n = 3. We assume that 0 C ffi 3 is a convex domain with the smooth 
boundary S. We assume that T is a connected and simply connected subset of S with the 
smooth boundary 7 . We will analyze T when x vanishes to a finite order k on 7 . With a slight 
abuse of notation, we arc-length parametrize 7 by the function 7 : R —>• M 3 (with 7 (s + L) = 7 (s), 
where L is the length of 7 ). 

Similarly to the case n = 2, we define 

Ay = {(x,£;x,Q : (x,£) £ V}. 


We denote by A the following canonical relation in (T*0 \ 0) x (T*P \ 0): 


A = {(x,t {x - 7 ( 5 )); y,r {y - 7 ( 5 ))) : \x - 7 (s)| = | y- 7 (s)|, 

(x — 7 (s), 7 / (s)) = (y — 7 (s), 7 / (s)) , for some s,t £ M such that r 7 ^ 0 }. 


We notice that (x, £; y,rj) £ A if and only if ( x , £) and (y, 7 ) are boundary elements corresponding 
to a common boundary point z £ 7 and they are obtained from the other by a rotation around 
the tangent line of 7 at z. 

The following result gives us a geometric description of the singularities of the Schwartz kernel 
/i of T. 


Proposition 3.1. We have 


WF{fi)' C Ay U A. 


The proof of Proposition 13.II is similar to that of Proposition [2T] (see also |Ngul5b Proposition 
3.1] and |FQ15| ). We skip it for the sake of brevity. Similarly to Proposition 12. 1 1 we obtain the 
following implications of Proposition 13.11 


a) T smoothens out all the visible singularities. 

b) T may reconstruct the visible singularities, and 

c) T may generate artifacts by rotating a boundary singularity around the tangent line of 7 
at the corresponding boundary point. 

The following result tells us the strength of the reconstructed singularities, explained in b), 
and artifacts, described in c): 

Theorem 3.2. The following statements hold: 

a) Microlocally on A \ A, we have y £ 7° (A) with the principal symbol 

o'ofoO = \ bcM + x(z-)], 

where z± is the intersection of the ray R±(x, f ) with S. 
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b) Microlocally on A \ A, we have /i£ / 4 2 (A). 

Similarly to Theorem 1 2. 31 we obtain the following implications of Theorem l3.2l (see also Lemma lA.181 
for the discussion on artifacts): 

a) If (x,f) £ WF(/) is a visible singularity, then (x ,£) is reconstructed with correct order. 

b) The artifacts are at least k order(s) smoother than the strongest generating singularity. 

Let us now proceed to prove Theorem 13.21 We will need the following two lemmas: 

Lemma 3.3. Let T : £'(Ll) —> be defined by 

= J J e i ^ z -\ 2 -\y-^ x a(x,y,z 1 ,X)f(y)dz J dydX, 

7 

where a £ S m ((Ll x fl x 7) x M). Then, T £ I m ~%( A). 

We first notice that the above formula J-(f) does not directly define an FIO, since the 
“phase” function (f> = (\x — z 7 | 2 — |y — z 7 | 2 )A involves an extra variable z 7 £ 7, which is neither 
a variable of / nor a phase variable. In |Ngul5b], where 7 is a line segment, the above result 
was proved by a change of variables. When 7 is a general curve, such change of variable seems 
to be complicated. We, instead, introduce a simple idea of lifting up the space. 

Proof. For the notational ease, let us denote X = Ll and Y = Ll x 7. Then, X and Y are smooth 
manifolds of dimensions nx = 3 and ny = 4, respectively. 

Let us define the operator C : £'(f 1 ) — > £'{fil x 7) by 0 

£(/) 0 , zfi) = f(x), (x, zfi) £ ft X 7. 

Then, C can be written in the following form 

£(f)( x ,Zy) = ^33 J J f(y) d V d C, (x,z 7 )£flx 7. 

R 3 R 3 

That is, C is an FIO of order (see (|28l) in Appendix IA. 310 

m c = (2 -n x ~ n Y )/4 = -i, 

with the canonical relation 


A c = {(z,%,£,0;x,£) : ( x ,£) £ T*fl\0j. 

4 This operator is just to lift up the space by one dimension. 
s Here, ng = 3 is the dimension of the phase variable £. 
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Let us define Jo : £'(T) —» D'(X) by the formula 

V r o{g){x) = J J e l tt x ~ z ^ 2 -\ y - z ^ 2)x a(x,y,z 1 ,\)g(y,z 1 )dz 1 dyd\. 

n 7 

Then, Jo is an FIO of order (see (f28l) in Appendix I A. 3 1) H 

mo = m + (2 n\ — nx — ny)/A = m — 

with the canonical relation 

A 0 = {( x , t ( x - z 7 ); y,z 7 ,r (y-z 7 ),r (y - x,j)) : \x - 7 (s)| = \y — for some r / 0}. 

We observe that T = Jo o C, and A = Ao o A£. Therefore, see [Hor7li pp. 96], 

T € J' m£+m o( A ) = J m_ i(A). 


□ 


Proof of Theorem 13.21 The proof of a) is similar to that of Theorem 12.31 a). We skip it for 
the sake of brevity. We now proceed to prove b). 

Let {x*, £*; y*, g*) E A \ A. That is x* / y* and there is z* = y(t*) E 7 such that 

\x* - z*\ = \y* - z*\, (x* — y*,r*) = 0. (14) 

Here, r* = is the unit tangent vector of 7 at z*. Let r| be the unit normal vector of 7 

which is tangent to S and points inward to T. Let d be the metric on S and O C S be a small 
(bounded) neighborhood of £* such that for each z E O there exists uniquely z^ E 7 such that 

d (z,Zry) = min{d(z, z') : z! E 7 }. 

That is, each z E O can be unique parametrized by (^ 7 , <5 = d(z 7 , 2 )). By narrowing down O, if 
necessary, we may assume 

z —>• (z 7 , 5) 

defines a smooth map from OflT to 7 x [0, <5o], for some 5q > 0, whose Jacobian | J(z)| is bounded 
from below. 

Let c E such that c(z) = 1 near z = z*. Let us write 

T:=T C + T', (15) 


where 

T c := Brc x VIZ, T := B v (1 - c) X V1Z. 

6 Here, n\ = 1 is the dimension of the phase variable A. 
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Let us denote by /x c , y' the Schwartz kernels of T c and T' , respectively. Then fi = y c + y!. 

An argument similar to the proof of Proposition 13.1 1 isee also proof of Proposition 12. ID shows 
that 

WF(/x'). 

It now remains to analyze y c . It can be written in the form: 

Me (x,y) = 2^3 / / e i(| — |2 -^— |2)A A 2 (z-x,v s ) c(z) X (z)dz d X. (16) 

k r 

By changing the variable z £ O n T —> (z 7 , s), we obtain 

<5 

y c (x,y) = ^3 y J J e i ^ x - z ^^ 2 - ly - z ^’ S) ^ x h(z 1 ,S)d6dz 1 X 2 dX. (17) 

R 7 0 

Here, 


h(zry,5) = (z - x,v z ) x(z)c{z)\J(z)\ 1 

satisfies h{z 1 ,8) = 0 for s > 5 and h(z^,$) vanishes to order k at s = 0. Moreover, notice that 
by choosing O small enough, we can assume that 8 as small as we wish. Let us show that there 
is a neighborhood Q of (x* , y* ) € H x H such that 

<5 

e i{\x-z(z^S)My-z{z-y,S)\ 2 ) A ^ dS = e i{\x-z^-\y-z^) A ^ ^ ^ A ) ? 

0 

where a(x, y, z 1 , X) £ S~ k (Q x 7 x R). (18) 

Indeed, consider the phase function of the left hand side 

A) = (|x - z( 2 : 7 , 5)| 2 - |y-z(z 7 ,5)| 2 )A. 

Taking the derivative with respect to s, we obtain 

4> s (x, y, z 7 ,8, X) = 2 (x - y, z s (^ 7 , 8)) X. 


Therefore, 


<j>s(x, y, z ly 0, A) = 2 (x - y, t 2 (z 7 )) A, 

where t 2 (z 1 ) is the unit vector tangent to S and normal to 7 (pointing inward to T). From (11411 
and the fact that H is convex, we easily observe that 


(x*-y*,r 2 (z*))^ 0 . 


( 19 ) 


Therefore, by choosing 8 > 0 small enough, we may assume that 

4> s (x, y, z, 0, A) / 0 for all (x, y , 27 , 8 , A) € Q x 7 x [0, <5] x (R \ 0), 
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where Q is a small neighborhood of (x,y) G x fb Taking integration by parts, we obtain for 
all (x, y) G Q: 


d 

f e i(|a:-a(« 7 ,<5)| 2 -|2/-2(2 7 ,<5)| 2 ) A ^ _ e i(\x-z~,\ 2 -\y-z~,\ 2 ) A h(z^, 0) 

J ifis (*£j V i Z'y , 0 , A) 

0 

<5 

+ [ e K|z-*(2 7 A)l 2 Hy-z(> 7 ,c5)| 2 ) x d s ( Hz'yifi) A ^ 

J vi<^(x,y,z 7 ,<f, A)/ 


o 

Continuing the integration by parts 

* fc 

i{\x-z{z~ t ,8)\ 2 -\y-z{z 1 ,5)\ 2 ) A ^ ^ _ e i(h- 2 7 | 2 -|y- 2 7 | 2 ) A 

Z=0 
<5 


0 


Hi(x, y, Zy,0, A) 


o 


+ J e i(l*-*(*7,«)l a -lw-*(*7,«)l 2 ) A 5 s i? fc+1 (x,y,^ 7j( 5) d<y. 


o 

Here, Hi is homogeneous of degree — l — 1 with respect to A. That is, 

s 

J e *(| a; -2(2 7 ,5)| 2 -| 2 y-2(2 7 ,5)| 2 )A^ ) ^ d( j = e i(|x-2 7 | 2 -| 2 /-2 7 | 2 )A a ( X)y) ^ ) ^ A )_ 


Here, 


i(x,y,z 11 X) = Hi(x,y, z 7 ,0, A) + y, z 7 , A) 


z=o 


and 


0 

R k (x,y, Zl , A) = J e*[d 3; - 2 ( 2 7 ’' 5 )l 2 -l y_2: ( z 7 ’' 5 )l 2 ) _ (l 3;_2 7 l 2_ l y_2 7 l 2 )l A d s H k+ i(x, y, z 7 , <5) d<5. 


Using standard integration by parts as above, one can obtain 

i4 € S -fc-1 ((Q x 7 ) x M). 


Moreover, from the definition of Hi we get Hi = 0 for all 0 < l < k — 1. Therefore, 
a(x, y, z 7 ,X) = H k (x, y, z 1 , 0, A) + R k (x, y, z 7 , A). 

This finishes the proof of (fl8l) . 

Now let us write: 

Mc(x, y) = ^J J e^ x -^ 2 -\y- z ^ A a(x, y, z 7 , A) A 2 cfe 7 d A. 

R 7 

Applying Lemma 13.31 we obtain /i c | gG I _fc_ 2 (A). □ 
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4 Numerical demonstrations 


In this section we investigate reconstructions from limited view circular mean data in a series of 
numerical experiments. Here, we focus on the two dimensional problem for illustration reasons. 
We will show that our theoretical results from Section [2] directly translate into practice and, in 
particular, that they can be used to significantly improve the reconstruction quality. 


Experimental setup. In all of our experiments we consider the phantom consisting of a disc 
with radius 0.3 centered at the origin, see Fig. [5] |(a)[ where the image size of the phantom 
was chosen to be 2048 x 2048. In our experiments, we will numerically illustrate the visible 
and invisible singularities, as well as boundary singularities for different angular range. We will 
also numerically investigate the difference between singly and doubly visible singularities. It is 
therefore useful to keep in mind that all the singularities of the original phantom are located on 
the circle of radius 0.3 centered at the origin, and the directions of all singularities are given by 
normal (orthogonal) directions to the circle at the location of the singularity. 

In what follows, we assume that the limited view data of this phantom are collected on a 
circular arc of the form 

T;, = {z(s) = (cos s, sins) : 0 < s < 6 }. 

In all of our experiments we computed the circular means of the phantom analytically and, 
to obtain limited view data, we sampled the data on T* (first experimetn) and (second 
experiment). In each experiment, we chose the number of angular samples n a as well as the 
number of radial samples n r independently of the angular range as n a = n r = 2048. Given 
this data we then implemented and applied the reconstruction formula Tof = Bxr b 'P'R.f hi 
Matlab (cf. ([3])) where xr 6 is the characteristic function of r^. In this situation, no artifact 
reduction was performed. To incorporate artifact reduction into the reconstruction formula, we 
also implemented the modified reconstruction formula Tf = BxVlZf where x was constructed 
to be smooth in the interior Int(F;,) of Tf, and, at the same time, to vanish to some order k > 0 
at the end points of T 5 (cf. (j4]) and subsequent discussion). More precisely, in our experiments 
we consider the following construction. Let 


H(s ) 


s(b — s ) 
s(b — s) + e’ 


and h(s) 


H(s) 


( 20 ) 


Then, h G C°°((0,b)) and 

i) 0 < h(s ) < 1 for 0 < s < b, h( |) = 1 , 

ii) h vanishes to order 1 at s = 0 , 6 . 

Let x be defined as x( z ( s )) = h(s), then T is smoothing of order one. The parameter e > 0 
controls how close the function x is to the constant function 1. The closer e is to 0, the closer 
h is to the constant function 1 on (0,6), see Fig. [4] |(a)| That is, the smaller e, the closer is 
the function x 1° the constant function 1 on T. We also consider smoothing of orders 2 and 3. 
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The corresponding smoothing functions are defined through h 2 (s ) and /i 3 (s), respectively. In 
general, we set 

X k (z(s)) = h k (s), k = 1,2,3,.... (21) 

The graphs of the functions h k for k = 1,2,3 with e = 0.2 are shown in Fig. [4]|(b)| 




(a) Smoothing function h for different (b) Smoothing functions h k of orders 

values of e k = 1, 2, 3 

Figure 4: Smoothing functions h k for different values of e and for different smoothing orders at 
the ends of the angular range, cf. (HQD-dUD. 


Experiment 1: In our first experiment we assume the data are collected on the acquisition 
surface (quarter of the unit circle) and consider reconstructions without as well as with 
artifact reduction (using the smoothing function (1201) with parameter b = ^). We numerically 
investigate the effect of changing the reconstruction parameters e as well as the effect of changing 
the smoothing order k at the ends of the angular range. The results of this experiment are 
presented in Fig. [MU 

Before examining the reconstructions, let us first note that the visible singularities at the 
considered angular range [0, are located on two open quarters of the circle: first (upper 
right) and third (lower left). Accordingly, the invisible singularities are located on the other 
two quarters of the circle, whereas the boundary singularities correspond to the boundary of the 
acquisition surface Tz and, therefore, are located at the four points (0,±0.3), (±0.3,0). These 
characterizations are due to the analysis and explanations that we presented in Section [2j To 
numerically verify those theoretical findings, let us now examine the reconstruction in Fig. [5j 
Our observations are as follows: 

a) All visible singularities are reconstructed sharply. They visually appear to be of the same 
order as the original singularities (jump from red to blue). 

b) The invisible singularities are smoothed away and, hence, not present in the reconstruction. 
This can be seen from the fact that there are no sharp boundaries (intensity jumps) along 
invisible directions. 

c) Added singularities (artifacts) are generated along four circles, each of them touches the 
disc (phantom) tangentially. Moreover, we observe that two circles are concentric and 
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(b) Reconstruction 



Figure 5: Reconstruction from limited view data, collected on T* (quarter of the unit circle), 
using the standard reconstruction operator 7o with no artifact reduction. The acquisition surface 


T| is illustrated by the green line in the phantom image 


centered on the x-axis, and the other two are concentric and centered on the y-axis. More 
precisely, the added artifacts are located on circles that are centered at the boundary 
points of the acquisition surface T* (which is illustrated by the green curve in Fig. [ Jl'a)[ ) 
and that are tangent to a singularity of the original phantom. That is, the artifacts are 
generated by the boundary singularities at (0.3,0), (—0.3,0), (0,0.3), and (0, —0.3). By 
further examining the artifacts in Fig. [ ^b)[ we also observe that the jumps along the 
added artifact circles are not as sharp as in the case of visible singularities. This indicates 
that the added artifacts are weaker than original (generating) singularities. In fact, our 
theoretical analysis shows that they are ^-order weaker. 

Summing up, this experiment shows that the above observations correspond to our theoretical 
findings stated in Propostion 12.11 and Theorem 12.31 

In the next step we investigate performance of the modified (artifact reduction) reconstruction 
operator T. To that end, we use the operator T with the cutoff function x k defined in (12U1) - (12TT) 
(cf. also Fig. [I]), and apply it to the limited view data. Note that the cutoff function x k is 
smooth in the interior of T* and vanishes to an order k at the end points of According 
to Theorem 12.31 the reconstructions obtained through T, will exhibit added artifacts that are 
k + \ orders smoother than the original singularities. Therefore, the degree of artifact reduction 
is linked to the order k and we expect the operator T to mitigate artifacts more when the bigger 
the order k is. In addition to that, we expect that the strength of artifacts is influenced by the 
parameter e (see (1201) and the definition of h k ). 

To investigate the practical performance, we computed a series of artifact reduction recon¬ 
structions by varying the parameters e and k. The results are shown in Fig. [6] and [7] First, 
we observe that in all reconstructions shown in Fig. [6] and [7] most of the visible singularities 
are well reconstructed. In Fig. [6l we have displayed some reconstructions using smoothing 
order k = 1 and varying the parameter e. Here we observe that for e = 0.05 almost no artifact 
reduction happens. This is due to the fact that, in the discretization regime, x changes very 
fast near the endpoints of Tjl and, hence, x behaves like a discontinuous function. The artifact 
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Figure 6: Reconstruction from limited view data, acquired on r* (quarter circle), using the 
modified reconstruction operator T (cf. (|4)) ) with the smoothing function x 1 defined in m- 
m- The figures illustrate the influence of the parameter e on artifact reduction. 


reduction gets clearer as we increase the value of e. Next, we consider the effect of varying the 
smoothing parameter k for a fixed e = 0.2. The corresponding reconstructions are shown in 
Fig. [U As expected, the artifacts get weaker (are better reduced) as the order increases. This 
is in accordance with our theoretical characterizations in Theorem 12.31 
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(a) k = 1, e = 0.2 
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(b) k = 2, e = 0.2 



Figure 7: Reconstruction from limited view data, acquired on T* (quarter circle), using the 
modified reconstruction operator 'T (cf. Q) with the smoothing function x k defined in (f20|) - 
(121)) . The figures illustrate the influence of the smoothing order k on artifact reduction for a 
fixed parameter e = 0.2. 


Experiment 2. Our second experiment follows the lead of our first experiment. Here, we only 

consider a larger angular range where we use limited view data collected on T 3tt (three quarters 

2 

of the unit circle). Again, we compute a series of reconstructions using the standard as well as 
the modified reconstruction operators, To and T, respectively. The results of this experiment 
are shown in Fig. 151-1101 
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Before we start, let us first remark that in this example all singularities of the phantom 
image are visible (they are located on the circle centered at the origin of radius 0.3), and the 
locations of all boundary singularities are the same as in Experiment 1, namely (±0.3,0) and 
(0,±0.3). In contrast to Experiment 1 where all of the visible singularities were singly visible, 
we now have both types of visible singularities, doubly and singly visible ones. Those on the 
first and third quarters are doubly visible, while those are on the second and fourth quarters are 
singly visible. See Section [2] for theoretical explanation. 


200 
400 
600 
800 
1000 
1200 
1400 
1600 
1800 
2000 

500 1000 1500 2000 500 1000 1500 2000 

(a) Original phantom (b) Reconstruction 



Figure 8: Reconstruction from limited view data, collected on T^r (three quarters of the unit 

2 

circle), using the standard reconstruction operator 7 q with no artifact reduction. The acquisition 


surface T 3 tt is illustrated by the green line in the phantom image (a) 
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(a) k = 1, e = 0.05 
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(b) k = 1, e = 0.2 
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(c) k = 1, e = 1 




Figure 9: Reconstruction from limited view data, acquired on (three quarters of the unit 

2 

circle), using the modified reconstruction operator T (cf. (JTJ) ) with the smoothing function x 1 
defined in (|20]) - ([2T|) . The figures illustrate the influence of the parameter e on artifact reduction. 


By examining the reconstructions using the standard reconstruction formula 7o (be. x — Xr) 
in Fig. [8] we easily observe that indeed all singularities of the phantom are reconstructed reliably. 
Here, the doubly visible singularities have more contrast than the singly visible ones. This is 
due to the fact that for each doubly visible singularity there are two positions on the acquisition 
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Figure 10: Reconstruction from limited view data, acquired on Tstt (three quarters of the unit 

2 

circle), using the modified reconstruction operator T (cf. J4])) with the smoothing function 
defined in (I20I) - (I21|) . The figures illustrate the influence of the smoothing order k on artifact 
reduction for a fixed parameter e = 0.2. 


arc r 3 tt from which the singularity is visible, whereas there is only one position on this arc for a 
2 

singly visible singularity. Mathematically, this is reflected by the different values of the principal 

symbol oo(a;,£) of the reconstructions operator T, cf. Theorem 12.31 where we can see that the 

principal symbol <7 o(x,£) = 1 if (x,£) is doubly visible and cr(x,£) = | if (x,£) is singly visible. 

In Fig. [9j we further observe that added artifacts are generated on circles that are centered at the 

boundary points of r^r and tangent to the boundary singularities. These artifacts, however, are 
2 

not as strong as the reconstructed singularities, which is again in accordance with our theoretical 
results, see Theorem 12.31 and the discussion below. We again studied the performance of artifact 
reduction by using the modified reconstruction operator 7~ with the smoothing function for 
b = (cf. (I20IMI21D ). The reconstruction results for varying parameters e and for varying 
smoothing orders are shown in Fig. [9] and Fig. [TOl respectively. Not surprisingly, we observe 
here the same behavior as in Experiment 1. 


Experiment 3. In our last experiment we investigate how the choice of the smooth cutoff 

function for T influences the artifact reduction. To that end, we consider the same limited view 

situation as in the Experiment 2 where the data are collected on T and define a new smoothing 

2 

function y which is equal to 1 in the interior of T 3^ and smoothly decreases to 0 in transition 

2 

regions of length e at the boundary of T 3^. To that end, we let 

2 


and 


h 0 (s) = -ys( 2e - s), 
e z 


0 < s < e 



' h 0 (s), 

1 , 

, ho(—s + 1), 


0 < s < e, 
e < s < 1 — e, 
1 - e < a < 1. 


( 22 ) 
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The parameter e > 0 again controls how close the function x is to the constant function 1. 

Moreover, this function h vanishes to order 1 at the endpoints of T^r. To obtain higher order 

2 

smoothness at the endpoints we again consider integer powers of h and set 

Xn ew(*(s)) = h new( s )- ( 23 ) 

A plot of the functions h k (s) is depicted in Fig. [11] for different values of e and k. The corre¬ 
sponding limited view reconstructions are presented in Fig. IT21 and [TUI 

The advantage of such a choice of the function x lies in the fact that it is exactly (not 
approximately) equal to 1 in the range 4ye < s < 4^(1 — e). Therefore, if e is very small, most 
of the visible singularities are reconstructed up to the factor 1 (if it is doubly visible) or ^ (if it 
is singly visible). However, the (theoretical) disadvantage of such a choice is due to singularities 
of x at the interior points z(^-e) and z(4f-(l — e)). According to our analysis in Section [2] this 
may lead to the generation of added artifacts (located on circles that are centered around these 
points). However, since x k is & order smoother at these inner points than at the endpoints, those 
new artifacts will be weaker than those rotating around the endpoints 2(0) and z(^-). Indeed, 
as can be seen in Fig. [T2l and [T3l the new added artifacts are too weak to be recognized in the 
reconstructions. 

Concerning the influence of parameters e and k, we arrive at similar observations as in Ex¬ 
periments 1 and 2. Comparing the reconstructions that were computed with different smoothing 
functions in Fig. fTOl and Fig. [13] we observe that the new smoothing function leads almost al¬ 
ways to a significantly better artifact reduction. For example, in Fig. dfsE the artifacts almost 
completely vanish and the phantom as well as the background are reconstructed very well. This 
examples shows that the choice of the smoothing function might influence artifact reduction 
performance significantly. 




(a) /iLw for e £ {0.1,0.25,0.4} 


(b) /in ew for £ = 0.25 and k = 1, 2, 3 


Figure 11: Alternative smoothing function h k ew for different values of e (length of smooth 
transition region at the boundary of the angular range) and for different smoothing orders k. cf. 

m-m- 


Conclusion. The above numerical experiments show that our theoretical results directly trans¬ 
late into practical observations. In particular, the proposed artifact reduction technique can lead 
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Figure 12: Reconstruction from limited view data, acquired on r^r (three quarters of the 

2 

unit circle), using the modified reconstruction operator T (cf. (pfj)) with the new smoothing 
function Xnew defined in (f22j) - (f23|) . The figures illustrate the influence of the length of the 
smooth transition region (parameter e) on artifact reduction for a fixed smoothing order 1. 


to a significant improvement of the reconstruction quality if the smoothing function is chosen 
appropriately. We will explore more experiments and report in-depth results in a future publi¬ 
cation. 
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A Appendix - Background in microlocal analysis 

Let P C M n be an open set, and T*P be the cotangent bundle of P. For simplicity, we can 
consider T*P as P x M n . We also denote 

T*P\0 = {(®,f) eT*P:^ 0}. 
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Figure 13: Reconstruction from limited view data, acquired on F^ (three quarters of the unit 

2 

circle), using the modified reconstruction operator T (cf. Q) with the new smoothing function 
Xnew defined in (f22|) - (f23|) . The figures illustrate the influence of the order of smoothing k on 
artifact reduction for a fixed length of the smooth transition region (e = 0.4). 


Let T>(Fl) = Cq°(FI) and V'(Q) be the standard spaces of test functions and distributions on 
fi. In this section, we briefly introduce some basic concepts in microlocal analysis, such as 
wave front set, pseudo-differential operators (4/DO), and Fourier Integral Operators (FIOs). 
Extensive presentations can be found in |Hor71l IHor831 ITreSOal ITVeSObj . The use of microlocal 
analysis in geometric integral transforms are pioneered in |Gui761 IGU89 I. iGU90bl Qui931. Its 
extensive uses in the studies of spherical mean transform can be found in many works, see 
|LQ00[ IXWAK091 QRS11 Ngul3 FQ15 Ngul5b| , just to name a few. 


A.l Wave Front Sets 

Definition A.l (Wave Front Set |Hdr71| ). Let f E D'{ Q) and (xo>£o) € T*fl \ 0. Then, f is 
microlocally smooth at (xo;£o) if there is a function <p E satisfying ip(x o) ^ 0 and an 

open cone V containing £o, such that J 7 ^pf) is rapidly decreasing in V. That is, for any N > 0 
there exists a constant Cn such that 

\Hvm)\ < C N ( 1 + \Z\)~ N , for all £ € V. 

The wavefront set of f, denoted by WF(f), is the complement of the set of all (a;o>£o) E T*fl 
where f is microlocally smooth. 


An element (x, £) E WF(/) indicates not only the location x but also the direction £ of 
a singularity of /. For example, if / is the characteristic function of an open set O <s with 
smooth boundary dO , then (x,£) E WF(/) if and only if x E dO and £ is perpendicular to 
the tangent plane of dO at x. Detailed discussion can be found in }Pet83| and, more briefly, in 

mm . 

Let T be a bounded linear operator from to T>'(Ll 2 ). The following rule provides an 

estimate of WF(/) in terms of WF (Tf), see |Hor71l Theorem 2.5.14]: 
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Theorem A.2. Let p be the Schwartz kernel ofT. Assume that WF(/j,) C (T*fl 2 \ 0 ) x (T*f2i\0) 
then 


WF(T /) C WF(n)' o WF(f). 


In the above theorem and elsewhere, WF(/j)' is the twisted canonical associated to WF(/r) 


WF(/i)' = {(x, £; y, -r,) : (x, £; y, rj) G WF( M )}, 


and 


WF( M )'oWF(/) := {(x,£) : (x^r?) G WF(/x)', for some (y, 77 ) € WF(/)}. 

The following theorem provides the product rule for wave front set, see [Hdr71 . Theorem 
2.5.10]: 

Theorem A.3. Let u,v be two distributions on f l. Then the product uv is well defined unless 
(x,£) G WF(u) and (x,—£) G WF(v) for some (x,£). Moreover, 

WF(uv) C {(x, £ + rf) : (x, £) G WF{u) or £ = 0, (x, rj) G WF(v) or 7 ] = 0}. 

The following theorem provides the composition rule for wave front sets, see jHor71l Theorem 
2.5.15]: 

Theorem A. 4. Let T\ and I 2 be linear transformations whose Schwartz kernels are ji\ G T>' (Ft \ x 
LI 2 ) and G V{Fl 2 x FI 3 ). We assume that WF(ni) C \ 0) x (¥*^2 \ 0) and WF(/j 2 ) C 

(T*n 2 \ 0) x (T*f2 3 \ 0). Then, the Schwartz kernel ji of T\ o T 2 satisfies: 

WF{n)' C WF(fi i)'o WF(fi 2 y. 

The following definition helps to quantify the strength of a singularity: 

Definition A.5 (Sobolev Wave Front Set [Pet83| ). Let f G Vifil) and (xo,£o) £ T*f2\0. Then 
f is in the space H s microlocally at (xo,£o) if there is a function <p G C{fi(Fl) satisfying <p(xo) 0 
and a function «(£) homogeneous of degree zero and smooth on 1 " \0 with u(£ 0 ) 0 , such that 

«(o^/)(e)€L 2 (K n ,( i + iei 2 r). 

The FP-wave front set of u, denoted by WF s (u), is the complement of the set of all (xo,£o) G 
T*fl where u is not microlocally in the space H s . 

One can use the Sobolev orders to compare the singularities (xi,£i) G WF(/i) and (x' 2 , £ 2 ) G 
WF(/ 2 ), where f\, f 2 are two distributions, not necessarily defined on the same set. For example, 
(xi, £ 1 ) is stronger than (X 2 , £ 2 ) if there is s such that (xi, £ 1 ) G WF s (/i) but (X 2 , £ 2 ) & WF s (/ 2 ). 

We also introduce the definition of conormal distribution (e.g., [Hdr711 lGIJ90al IFLU031 
IEswl2j l: 
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Definition A.6. Assume that S C ft is a smooth surface of co-dimension k. Let h £ C°°(fl,lR k ) 
be a defining function for S with rank(dh ) = k on S. The class I r (S ) consists of the distributions 
which locally can be written down as a finite sum of oscillatory integrals of the form 

* ) = f **»*,*)», 

R k 


where a £ S r (fl x R fc ). 

In the above definition and elsewhere in this article, we use the following definition of a 
symbol: 

Definition A.7 f |H6r71j ). Let Q C R n be an open set. The space S m (Ll x R^) consists of all 
functions a £ C°°(Ll x (l w \ 0)) such that for any multi-indices a, fi and K <e f l, there is a 
positive constant C a ^ t K such that 

\&Zd$a{x,£)\ < C aAK ( 1 + |e|) m “ H , for all (x,0 £ K x (1^ \ 0). (24) 

The elements of S m (Ll x 1^) are called symbols of order m. 

We note that if u £ I r (S ), then WF(u) C N*S (see, e.g, [Hor71| ). where N*S is the conormal 
bundle of S. 

Definition A.8. Let f £ D'(R 2 ). We say that (xo>£o) € WF(f) is a conormal singularity of 
order r to the surface S if there is u £ L r (S) such that 

(®o,£o) ^ WF(f-u). 

One can use the order r to compare two conormal singularities (x,£) £ WF(/i) (along the 
surface Si) and ( y,rj ) £ WF(/ 2 ) (along the surface £ 2 ), where /i ,/2 are two distributions on hi. 
For example, (x,£) is weaker than ( y,rj), if there is r £ M such that (x,£) is of order r while 
(y, rj) is not. 

A.2 Pseudo-Differential Operators ('I'DOs) 

Given a £ S m (Ll x M n ), the operator T : Cq°(LI) -> C°°(P) defined by the oscillatory integral 

77 ( x) = J j a(x, 0 f(y) df dy, (25) 

Q R n 

is called a pseudo-differential operator ('FDO) of order m with the (full) symbol a(x,f). 

Since the above integral may not converge in the classical sense, the expression in (I25p needs 
to be properly defined, see, e.g., (Hor71i Proposition 1.1.2 ]. Given this proper definition, T 
extends continuously to £'(Ll) —> T>'(Ll). In particular, it can be shown that a pseudo-differential 
operator T does not generate new singularities. That is, |Hdr711 Page 131] 

WF(T/) C WF(/). (26) 
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Moreover, if / is in the space H s microlocally at (xo>£o) then Tf is in space H s m microlocally 
at the same element (xo,£o)) see I *cl 83 . Mi:, 

We will denote p £ I m ( A) if p is the Schwartz kernel of a pseudo-differential operator of 
order m. Let us define a technical term that is used in the statement of Theorem 12.31 ah 

Definition A.9. Let A C A be a conic set that is open in the induced topology of A. We say 
that near A, p is microlocally in the space I m ( A) with the symbol a if the following 
holds: for each element (x*, £*; x*, £*) £ A there exist p* £ I m (A) such that 

(z*,C;z*,n 0 WF(p-p*y, 

and the symbol of p* is equal to cr(x,£,) in a conic neighborhood if (,t*,£*). 

The following result was proved in |Ngul5a| , which is used to explain our result in Sec¬ 
tions [2] & [3) 

Lemma A.10. Let T : £'(Ll) —> T>'(Ll) be a linear operator whose Schwartz kernel p £ V(Q x Ll) 
satisfies WF(p)' C (T*fl\0) x (T*fl\0). Assume that p is microlocally in I m { A) near A with the 
symbol <r(x,£). Let ( x *,£*) £ T*fl\0 such that (x*, x*, £*) £ A and in a conic neighborhood 

of(x*,n 

k(*,oi>c(i + i£ir, ki>i. 

Assume further that 

{(x*,e-,y,v) e WF(p)' : (y,ri) € WF(f)} C A. (27) 


Then, for any s £ R, 

(x*,C) G ^(Z) if and only if (z*,£*) € WF s _ m {Tf). 

A.3 Fourier Integral Operators (FIOs) 

In this section, we introduce some special Fourier distributions which are needed in this article. 
The reader is referred to, e.g., [ Hor7L [Tre80bl Dui 1 I for the general theory of the topic. Let X 
and Y be two manifolds of dimension nx and ny, respectively, and A be a homogeneous canonical 
relation in (T*A \ 0) x (T*y \ 0). Then, there is an associated class of Fourier distributions of 
order m, denoted by I m ( A). Each element of I m ( A), called a Fourier integral distribution of 
order m, is a distribution p £ V(Fl x fi) such that it can be locally written down in the form 

/a(x, y) = J e 1 ^’^ a(x, y, A) dX. 

R N 

Here, 0 

a(x, y, A) £ s m+ ^ nx+nY - 2N)/A {X x Y x R N ), (28) 

and 0 is a phase function associated to A. That is, <f> = 4>(x,y, A) € C°°(X x Y x (R^ \ 0)) 
satisfies 

7 The order m + (nx + ny — 2iV)/4 of a specified here is due to, e.g., |Hor71l pp. 115]. 
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1) (j) is homogeneous of degree 1 in A, 

2 ) 4 > x 7 ^ 0 and 4 > y ^ 0 on the set 

C = {(x, y,X) e X xY x (R n \ 0) : d x cj) = 0}, 

and 

{(x, d x <t>\ y, -d y 4>) : (x, y, A) E C} C A. 

The following result gives us a rule for the wave front set of a Fourier integral distribution (see, 
e.g., |Hor71 . Theorem 3.2.6]): 

Theorem A.11. Let y E I m { A) then 


WF(y)' C A. 

The linear operator T : T>(Y) —> T>'(X) whose Schwartz kernel is y is a Fourier integral operator 
(FIO) of order m. With a slight abuse of notation, we also write T E J m (A). 

The following technical term is used in the statement of Theorem 12.31 b): 

Definition A.12. Let A C A be an open conic set in the induced topology of A. We say that 
near A, y is microlocally in the space L m ( A) if the following holds: for each element 
(x*,^*-,y*,ri*) E A there exists y * € I m ( A) such that 

(x*,£*;yW) 2 WF{y — y*)'. 

The following result (see |Hor71l Theorem 4.3.2]) is used to analyze the mapping properties 
of the FIOs discussed below: 

Theorem A. 13. Let O, C and A C (T*f2\0) x (T*fl\0) be a homogeneous canonical relation 
such that both of its left and right projections on have surjective differentials. Assume that 
the differentials of the left and right projections A —y TT2 have rank at least l + n. Then, every 

s m n ~ l 

T E I m ( A) maps continuously from Hf omp (Ll) to H loc 2 (Ll). 

In the discussion below, we introduce two classes of Fourier distributions whose canonical 
relation is defined by rotations around a point or around tangent lines of a smooth curve, 
respectively. 

A.3.1 Fourier distributions associated to a point 

Let us now introduce the class of Fourier distributions whose canonical relation is defined by 
the rotations around a point. For this type of distribution, we only consider Ll C R 2 . 

Let xo E ffi 2 such that xq ^ LI. We define the following homogeneous canonical relation in 
(T*\0) x (T*fi\0) 

A^o = {(^ T i x ~ x o); y, r{y - y 0 )) € (T*f) \ 0) x (TT2 \ 0) : |x — x 0 | = | y- x 0 |}. 
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That is, A Xo is defined by rotating (y, r) = r(y — xq)) € T*0 \ 0, which pass through xo, around 
x 0 . 

In Section [2l we work with the following explicit form of an element y € I m (A xo ): 

H{x,y) = J e i(l*-*o| 2 -|w-*o| 2 )A a (x,y,X)dX, 

R 

where a € S m+ i (Q x Q x I). 

The following result is a direct consequence of Theorem I A. 131 

Lemma A.14. Let T € I m { A xo ). Then, T maps continuously from Hf omp {Lt) to 

Proof. We only need to apply Theorem IA. 131 with 1 = 1. □ 

The following result is used to analyze the strength of artifacts in Section El Its proof is 
almost exactly the same as that of |Ngul5a[ Corollary 2.15]. We skip it for the sake of brevity. 

Lemma A.15. Let T : £'(11) —> T>'(Lt) be a linear operator whose Schwartz kernel y € T>'{LL x fl) 
satisfies WF(y) C (T*f2 \ 0) x (T*fl \ 0). Assume that y is microlocally in I m ( A) near an open 
conic set A C A. Let (x*,f*) € WF(Tf) D7rL(A) such that 

{{x*, y, rj) € WF(y,)' U A : (y,rj) € WF(f)} is a compact subset of A. 

If (a;*,£*) € WF S {Tf), then there is ( y*,rj*) € T*fi n A such that 

e A and (y*, V *) € WF s+m+ i(f). 


The following result is useful to analyze the artifacts when the original singularities are 
conormal. Its proof is almost exactly the same as that of |Ngul5a Theorem 2.16]. We skip it 
for the sake of brevity. 

Lemma A.16. Suppose that all the assumptions in Lemma \A.15 1 hold. Assume further that: 

1) There are at most finitely many y* G Ll such that 

(x*,t.*',y*i y* = r(y* — x 0 )) € A and (y*, V *) € WF(f). 


2) Each such (y*,rj*) is a conormal singularity of order r along a curve C whose contact order 
with S(rco, \x — xo] = |y — xo|) exactly 10 

Then, (x*,f*) is a conormal singularity of order at most m + r along the circle §(xq, \x* — xq])- 


8 We note here that the contact order is always at least 1, since both curves are perpendicular to g* at y*. 
Therefore, the condition on the contact order is quite generic. 
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A.3.2 Fourier distributions associated to a smooth curve 

Let us consider C R 3 . We introduce a class of Fourier distributions, whose canonical rela¬ 
tion is defined by the rotations around tangent lines of a smooth curve. This class of Fourier 
distributions appears in the statement and proof of Theorem 13.21 b). Let 7 be a closed smooth 
curve in R 3 parametrized by the parameter s. Assume that 7 fl fl = 0. We define the following 
homogeneous canonical relation in (T*Q \ 0) x (T*17 \ 0) 


A 7 = {(x, r (x — 7 ( 5 )); y,T (y — 7 (s))) : |x- 7 (s)| = |y~ 7 (s)|, 

(x - 7(s),t'(s)) = (y -'y{s),V(s)) ,x £ n,y € Q,,s € K,0 / t € R}. 


That is, A 7 is defined by rotating an element (y,r(y — 7 ( 5 ))), that passes through z = 7 (s), 
around the the tangent line of 7 at In Section [31 we make use of this class I m ( A 7 ). We state 
here some needed basic facts of this class. 

The following property is a direct consequence of Theorem IA. 131 (for 1 = 1): 


Lemma A.17. Assume that T £ J m (A). 


H. 


loc 


m 


Then, F is a continuous map from Hf om ( fl) —>• 


We note that: 


7 Tl(A) = ttr(A) = {(x, r (x — 7 ( 5 ))) :xGO,sG1,0/tG M}. 

The following result is a microlocal version of the above result, which is used in Section [3] to 
analyze the strength of artifacts. Its proof is almost exactly the same as that of | |Ngul5a 
Corollary 2.15]. We skip it for the sake of brevity. 

Lemma A.18. Let T : S'(fl) —> 'D'(fi) be a linear operator whose Schwartz kernel y € T>'(Q x Q) 
satisfies WF(y) C (T*fl \ 0) x (T *£l \ 0). Assume that y is microlocally in I m { A) near an open 
conic set A C A. Let (x*,£*) € WF(Tf) n 7 t^(A) such that 

{(x*, y, if) £ WF(y)' U A : (y, rf) £ WF(f)} is a compact subset of A. 

If (x*,£*) £ WF S (Tf), then there is ( y*,rj*) £ 7Tr(A) such that 

(x*,C;y*,y*) e A and ( y *, v *) e WF s+m+ i(f). 
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